Analytic, Group-Theoretic Wave Functions for Confined TV-Body Quantum Systems 



W.B. Laing 1 , M. Dunn 1 , J.G. Loeser 2 , and D.K. Watson 1 
1 University of Oklahoma, Department of Physics and Astronomy, Norman, OK 73019 
2 Oregon State University, Department of Chemistry, Corvallis, OR 97331 
(Dated: February 2, 2008) 

Systems involving iV-identical interacting particles under quantum confinement appear in many 
areas of physics, including chemical, condensed matter, and atomic physics. We discuss a beyond- 
mean-field perturbation method that is applicable to weakly, intermediate and strongly-interacting 
systems. Group theory is used to derive an analytic beyond-mean-field correlated wave function at 
zeroth order for a system under spherical confinement. We derive the corresponding zeroth-order 
analytic density profile and apply it to the example of a Bose-Einstein condensate. 
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INTRODUCTION 

During the last two decades, novel ./V-body quantum 
systems have been created using techniques to confine 
and manipulate atoms, ions, and electrons. These sys- 
tems, of both fundamental and technological interest, in- 
clude condensed atomic Bose gases, atoms confined in op- 
tical lattices, quantum dots, and ultracold fermion gases. 
Tuning the external fields of these environments provides 
unique opportunities to study many-body effects over a 
range of interaction strengths. 

Mean-field treatments, such as the Hartree-Fock 
method in atomic physics and the Gross-Pitaevskii 
method in condensed matter, do not include correlation 
effects and fail to describe systems with tight confine- 
ment or strong interaction. These systems, which have 
hundreds to millions of particles, present serious chal- 
lenges for existing iV-body methods, many of which were 
developed with small systems in mind. 

The methodology described below uses dimensional 
perturbation theory (DPT) [l| , which has been previously 
applied to mostly small- AT systems in the form of high- 
order, largely numeric calculations. In this letter, we 
use group theoretic techniques to develop an analytic ap- 
proach which fully exploits the symmetry of the zeroth- 
order problem. This method avoids heavy numerical 
computation, and iV enters into the theory as a param- 
eter. This crucial simplification allows results for any N 
to be obtained in a single calculation. This method also 
directly accounts for each two-body interaction, rather 
than using an average interaction. Even lowest-order re- 
sults include beyond-mean-field effects. Thus, in contrast 
to the low-density expansion methods pioneered by Lee, 
Huang and Yang in the 195Os0, this method is appropri- 
ate for the study of both weakly and strongly interacting 
systems and the transition between them. This general 
formalism offers a systematic approach to the study of 
correlation in atomic and molecular, condensed-matter, 
chemical, and nuclear systems. 

Nearly all past work using DPT has focused on en- 
ergies with little attention given to wave functions. In 



this paper we derive an analytic, correlated lowest-order 
S-w&ve wave function for N identical particles in a spher- 
ical confining potential. The lowest-order wave function 
yields important information such as the nature of ex- 
citations and expectation values of physical observables. 
This result can be systematically improved by going to 
higher order. 

TOOLBOX 

The tools used to describe large-iV correlated wave 
functions are carefully chosen to maximize the use of 
symmetry and minimize the dependence on numerical 
computation. We handle the massive number of inter- 
actions for N large ( ~ N 2 /2 two-body interactions) by 
bringing together three theoretical methods. 

The first, DPT0, is chosen because its zeroth- 
order equation, which is obtained for large D, yields a 
maximally-symmetric configuration for N identical par- 
ticles. Higher orders yield insight into fundamental mo- 
tions as well as a framework for successive approxima- 
tions. The second method is the FG method of Wil- 
son, Decius, and Cross0. This seminal method has long 
been used in quantum chemistry to study vibrations of 
polyatomic molecules. It directly relates the structure 
of the Schrodinger equation to the coordinate set which 
describes the normal modes of the system. The third 
method, the use of group theoretic techniques 0,01, takes 
full advantage of the symmetry at zeroth-order. 

Dimensional Perturbation Theory 

For TV-body systems in large dimensions, the DPT 
wave function is localized about a symmetric structure in 
which each particle is equidistant and equiangular from 
every other particle. The Jacobian-weig htcd 5J DPT 
wave function is harmonic and corresponds to oscillations 
about this structure. Notwithstanding its relatively sim- 
ple form, the large-dimension, zeroth-order wave function 
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includes beyond-mean-field effects. 

The Schrodinger equation for the zeroth-order 
Jacobian-weighted wave function has the formQ 

{^\d r T G d r + ^y' T Fy' + v^j $ (y') = e $ (y'), 

(1) 

where y' is the displacement coordinate vector formed 
from dimensionally-scaled internal displacement coordi- 
nates f'i (1 < i < N) and % (1 < i < j < A). Coordi- 
nates r' and 7', are related to the Z?-dimensional scalar 
radii of the A particles from the center of the confining 
potential and the cosines 7^ of the A(A — l)/2 angles 
between the radial vectors by 

n = k(D) (foe + S^fty , and 7ij = 7^ + tfVa^. . 

(2) 

The expansion parameter is S = l/D, and k(D) 
is quadratic in D with the particular form chosen to 
simplify the resulting equations^. The quantities foo 
and 7^ define the symmetric large-dimension structure, 
which depends on the nature and strength of the interpar- 
ticle interaction and confining potential. The matrices G 
and F are constants derived from the Hamiltonian. 



The F G Method 

The FG method^] is used to obtain normal-mode co- 
ordinates and frequencies from the eigenvalue problem 
in Eq. The 6 th normal- mode coordinate q' b may be 
written as 

[q'] b = b T y', where F G b = A b b , 
b T Gb = l, \ b =wl 

and tub is the normal-mode frequency. Equation @ still 
represents a formidable eigenvalue problem unless A is 
quite small, since there are P = A(A + l)/2 normal 
coordinates and up to P distinct frequencies. 

Group Theory: the Sn Symmetry 

The full S-wave Hamiltonian is invariant under parti- 
cle interchange. This fact defines an Sn symmetry un- 
der which the system is invariant. As defined in Eq. J2J, 
the large-dimension structure is a completely symmet- 
ric configuration so Eq. Q is also invariant under the 
group Sn. This Sn symmetry .4] brings about a remark- 
able reduction from P possible distinct frequencies to five 
actual distinct frequencies and greatly simplifies the de- 
termination of the normal coordinates through the use of 
symmetry coordinates. [3| 



As a prelude to the above, we note that the Sn invari- 
ance of Eq. Q means that the F, G and FG matrices of 
Eq. J2J arc invariant under Sn, which implies that the 
eigenvectors b and normal modes transform under irre- 
ducible representations (irreps.) of Sn- Using the theory 
of group characters^, the coordinates f ■ are reduced to 
one 1-dimensional and one (N — l)-dimcnsional irrep., 
labelled [N] and [N - 1,1] respectively 0. The are 
reducible to one 1-dimensional, one (N — l)-dimensional, 
and one N(N — 3)/2-dimensional irrep., labelled [N], 
[N — 1, 1] , and [N — 2, 2] respectively. Since the normal 
modes transform under irreps. of Sn and are composed 
of linear combinations of elements of vectors f ' and 7' , 
there will be two 1-dimensional, two (N— l)-dimensional, 
and one entirely angular N(N — 3)/2-dimensional irreps. 
labelled by the partitions [N], [N - 1, 1] and [N - 2,2] 
respectively. All normal modes that transform together 
under the same irrep. have the same frequency, so rather 
than P distinct frequencies there are only fivel 

NORMAL MODES - THE PROGRAM 

We determine the normal coordinates and distinct fre- 
quencies in a three-step process @: 

First, we define sets of primitive irreducible coordi- 
nates that have the simplest possible functional form sub- 
ject to the requirement that they transform under par- 
ticular non-orthogonal irreps. of Sn- We define two sets 
of linear combinations of elements of the r' vector which 
transform under non-orthogonal [N] and [AT — 1, 1] irreps. 
of Sn ■ We then derive two sets of linear combinations of 
elements of the 7' vector which transform under exactly 
these same two irreps. of Sn- Finally we define a set 
of linear combinations of elements of 7 which transform 
under a particular non-orthogonal [A — 2, 2] irrep. of Sn- 

Second, we use linear combinations within each set of 
primitive irreducible coordinates to determine symmetry 
coordinates that are defined to transform under orthog- 
onal irreps. of Sn ■ Care is taken to ensure that this 
transformation to the symmetry coordinates preserves 
the identity of equivalent representations in the r' and 
7' sectors. We choose one of the symmetry coordinates 
to be a single primitive irreducible coordinate, the sim- 
plest functional form possible that transforms irreducibly 
under Sn ■ The next symmetry coordinate is chosen to be 
composed of two primitive irreducible coordinates, and so 
on. Thus the complexity of the symmetry coordinates is 
minimized, building up slowly as symmetry coordinates 
are added. 

Third, the FG matrix, which was originally expressed 
in the r' and 7 basis, is now expressed in symmetry co- 
ordinates. This results in a stunning simplification. The 
P x P eigenvalue equation of Eq. ® is reduced to one 
2x2 eigenvalue equation for the [A] sector, A — 1 iden- 
tical 2x2 eigenvalue equations for the [A — 1,1] sector, 
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and N(N — 3)/2 identical lxl eigenvalue equations for 
the [N - 2,2] sector. For the [N] and [N - 1,1] sec- 
tors, the 2x2 structure allows for mixing of the f' and 
symmetry coordinates in the normal coordinates (see 
Eq. © below). The lxl structure of the equations in 
the [N — 2,2] sector reflects the absence of r' symmetry 
coordinates in this sector, i.e. the [N— 2, 2] normal modes 



are entirely angular. 

The Symmetry Coordinates. 



Using steps one and two, we derive the symmetry co- 
ordinates 



J 



l[N] 



1 N 



[N] 



N l-l 



k=l 

1 



where 1 < i < N — 1 , and 
1 



i N 1 T i N 



U=2 fc=l k=l l=k+l 



U-l 



Z=i+2 



7//' 



t+1 ' 



where 1 < i < j — 3 and i + 3 < j < N . 



j'=i+2 



(4) 



The Normal Coordinates. 



The Wave Function. 



In the third step, Eq. © is expressed in the symmetry 
coordinate basis of Eq. J3J and reduces to three eigensys- 
tem equations ( a = [N], [N — 1], [N — 2, 2] ): 



,FG„a 



X a c a , and [c a Y cr^c" = 1 



(5) 



The wave function in Eq. (JTJ is the product of P 
harmonic-oscillator wave functions: 



$o(y') 



n n 

= {0±, 1±, 2} £=1 



(7) 



The reduced FG matrices er^ and cr^_ l 1 j are 2x2 
matrices, while o"^L 2 2 ] is a 1 x 1 matrix. The same 
is true for the diagonal reduced G matrices, ■ The 
elements of <t^ g and cr^ are known analytic functions [y]. 

There are five solutions to Eq. 10 denoted ± = 



{A 



[JV]> 



} 



= {A 



[JV-1,1]' 



} and 2 = 



{AfAr-2,21; c^ -2 ' 2 !} . The two-element c± vectors for the 
a = [N] and [./V — 1, 1] sectors determine the angular- 
radial mixing of the symmetry coordinates in a normal 
coordinate of a particular a . Hence 



[q'] 6 = [c£V [s« ] £ + [dZfr [s«,]« . 



(6) 



The normal coordinate label b is replaced by the labels a, 
£ and ± on the rhs of Eq. ©. For the [N-2, 2] sector, the 
symmetry coordinates are also the normal coordinates up 
to a normalization constant, [c^ -2 ' 2 !]^' . 



where <b n „ \\fUT u 



[q /Al ]^j is a one-dimensional harmonic- 
oscillator wave function of frequency ui^, and n^, is the 
oscillator quantum number, < n AI? < oo, which counts 
the number of quanta in each normal mode. The quantity 
fi labels the manifold of normal modes with the same 
frequency w M while — 1 , N — 1 or N(N — 3)/2 for 



— ± , 1 ± or 2 respectively. 



UNMASKING THE WAVE FUNCTION: THE 
DENSITY PROFILE 



In the case of macroscopic quantum-confined systems, 
such as a BEC, the wave function is made manifest in 
the experimentally accessible density profile. 

The large-dimension Jacobian-weighted ground-state 
density profile, No(r) , is an analytic function: 
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FIG. 1: Number density per atom versus radial distance of a 
spherically confined BEC of 100 87 Rb atoms with a = 1 000 
a.u. and uj^o = 2n x 77.87 Hz. The solid line is the analytic 
DPT density and the dashed line is the MGP density. 

47rTV (r) = 47r r 2 p (r) 

N pOO POO 

=E/ ••• / sir-nwow n n%'i« 

. =1 J_oo J-oo M =0±,l±,2 £=1 

= JvJ-exp[-ii(r-K(D)f 00 ) 2 ], (9) 

V 7T 

where i? depends on TV, on the frequencies <D ± and W]± , 
and on the angular-radial mixing coefficients [c± ]f' , 

[cT ] h> , [cf- 1 ' 11 ],' -d ■ 



The analytic DPT density profile of Eq. (9) is a Gaus- 
sian, symmetric about K(D)f 00 . For a or TV sufficiently 
large, the physical density profile develops an asymme- 
try which motivates future work on the next-order DPT 
wavefunction. 

SUMMARY AND CONCLUSIONS 

This letter discusses an analytic approach to the study 
of quantum-confined TV-body systems. Unlike mean-field 
methods, this approach directly accounts for each inter- 
particle interaction, even at lowest-order, and so is appli- 
cable to high-density, strongly-interacting systems. By 
taking advantage of three powerful methods, DPT, the 
FG method and the group theory of the Sn group, we 
avoid heavy numerical computation and offer a system- 
atic approach to correlation that is not limited to small- TV 
systems. While most prior work with DPT has focused 
on energies, this letter derives the lowest-order, analytic, 
correlated wave function for a spherically-confined, TV- 
particle system and the corresponding analytic density 
profile. 
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An Example: the BEC Density Profile 

We perform a proof-of-concept test on our zeroth-order 
density profile by considering a T = OK BEC in an 
isotropic, harmonic trap. We take the interatomic po- 
tential to be a hard sphere of radius a = 0.0433a; lo (ten 
times the natural scattering length of 87 Rb) and TV = 100 
(a ten-fold increase in TV beyond existing benchmark dif- 
fusion Monte Carlo (DMC) density results @)- At this TV 
and a the modified Gross-Pitaevskii method (MGP) 9] is 
still valid for comparison. To implement the perturbation 
theory with k(D) = D 2 a,h , we allow this hard-sphere 



potential to soften away from D 
teratomic potential 



3 by defining the in- 



(1-3/79) 



1 — tanh 



c a (ry - \/2a) 
V2(l - 3/D) 



.... 
where a = D 2 a, and the parameters Vq and cq specify the 

height and width. We exploit this freedom in the shape 

and position of the interatomic potential at large D to 

ensure that the lowest-orders energy and wave function, 

^o(y')i are as accurate as possible. The parameters V 

and c~o are determined by a least-squares fit to benchmark 

DMC energies (TV = 5 — 100) and densities (TV = 3 and 

10) 8J. In Figure n° ur analytic density profile compares 

well to the MGP result. 
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